Consignes

Complétez ce document en remplissant les chunks vides pour écrire le code qui vous a permis de répondre à la question. Les réponses attendant un résultat chiffré ou une explication devront être insérés entre le balises html code. Par exemple pour répondre à la question suivante :

La bioinfo c'est : <code>MERVEILLEUX</code>.

N’hésitez pas à commenter votre code, enrichier le rapport en y insérant des résultats ou des graphiques/images pour expliquer votre démarche. N’oubliez pas les bonnes pratiques pour une recherche reproductible ! Nous souhaitons à minima que l’analyse soit reproductible sur le cluster de l’IFB.

Introduction

Vous allez travailler sur des données de reséquençage d’un génome bactérien : Bacillus subtilis. Les données sont issues de cet article :

Analyses

Organisation de votre espace de travail

# Je travaille dans mon répertoire de travail qui est associé au projet "dubii2021" sur le cluster de l'IFB.
cd /shared/projects/dubii2021/glelandais/modules-4-5-evaluation/.
# Création de répertoires pour organiser les fichiers.
mkdir raw-data
mkdir quality-controls
mkdir clean-data

Téléchargement des données brutes

Récupérez les fichiers FASTQ issus du run SRR10390685 grâce à l’outil sra-tools [1]

cd raw-data
# Comme je vais utiliser les outils en dehors d'un script,
# il me semble que cette commande est importante.
salloc --cpus-per-task=10 --mem=1G
# Chargement du module
module load sra-tools
# Lancement du téléchargement
srun fasterq-dump -p --split-files SRR10390685
# Deux fichiers FASTQ ont été obtenus

Combien de reads sont présents dans les fichiers R1 et R2 ?

# Le nombre de reads est égal au nombre de lignes divisé 
# par 4 (un read comporte 4 lignes)
grep -c "^@" SRR10390685_1.fastq
grep -c "^@" SRR10390685_2.fastq

Les fichiers FASTQ contiennent 7066055 reads.

Téléchargez le génome de référence de la souche ASM904v1 de Bacillus subtilis disponible à cette adresse

wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/009/045/GCF_000009045.1_ASM904v1/GCF_000009045.1_ASM904v1_genomic.fna.gz

Quelle est la taille de ce génome ?

# Décompression du féichier génome
gunzip GCF_000009045.1_ASM904v1_genomic.fna.gz
# Décompte du nombre de base
cat GCF_000009045.1_ASM904v1_genomic.fna | grep -v "NC_" | tr --delete "\n" | wc -c

# Note : J'ai trouvé de l'aide ici :
# https://dridk.me/genome_chiffre_1.html

La taille de ce génome est de 4215606 paires de bases.

Téléchargez l’annotation de la souche ASM904v1 de Bacillus subtilis disponible à cette adresse

# Téléchargement
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/009/045/GCF_000009045.1_ASM904v1/GCF_000009045.1_ASM904v1_genomic.gff.gz
# Décompression du fichier
gunzip GCF_000009045.1_ASM904v1_genomic.gff.gz

Combien de gènes sont connus pour ce génome ?

# Récupération de la troisième colonne du fichier et
# dénombrement des occurences de gene
cut -f 3 GCF_000009045.1_ASM904v1_genomic.gff | grep -c gene

4536 gènes sont recensés dans le fichier d’annotation.

Contrôle qualité

Lancez l’outil fastqc [2] dédié à l’analyse de la qualité des bases issues d’un séquençage haut-débit

# Changement de répertoire de travail
cd ../quality-controls
# Chargement du module fastqc
module load fastqc
# Lancement des analyses qualités
srun fastqc ../raw-data/SRR10390685_1.fastq -o .
srun fastqc ../raw-data/SRR10390685_2.fastq -o .
# Rapport MultiQC
module load multiqc
srun multiqc -d . -o .

La qualité des bases vous paraît-elle satisfaisante ? Pourquoi ?

  • Oui
  • Non

car les valeurs de score qualité sont élevées (> Q30) comme le montre le graphique “Per base sequence quality”

Lien vers le rapport MulitQC

Est-ce que les reads déposés ont subi une étape de nettoyage avant d’être déposés ? Pourquoi ?

  • Oui
  • Non

car tous les reads n’ont pas la même longueur

Quelle est la profondeur de séquençage (calculée par rapport à la taille du génome de référence) ?

# Utilisation de la formule du cours :
# N*L/G ; avec N = Nombre de lectures, L = Longueur des lectures et G = Taille du génome 
# ((7066055 * 2) * 151) / 4215606 = 506.2021  

La profondeur de séquençage est de : 500 X.

Nettoyage des reads

Vous voulez maintenant nettoyer un peu vos lectures. Choisissez les paramètres de fastp [3] qui vous semblent adéquats et justifiez-les.

# Changement du répertoire de travail
cd ../clean-data
# Chargement du logiciel fastp
module load fastp
# Lancement de fastp
srun fastp --in1 ../raw-data/SRR10390685_1.fastq --in2 ../raw-data/SRR10390685_2.fastq --out1 ./SRR10390685_1_clean.fastq --out2 ./SRR10390685_2_clean.fastq --html ./fastp.html --cut_mean_quality 30 --cut_window_size 8 --length_required 100 --cut_tail --json ./fastp.json

Les paramètres suivants ont été choisis :

Parametre Valeur Explication

Ces paramètres ont permis de conserver reads pairés, soit une perte de % des reads bruts.

Alignement des reads sur le génome de référence

Maintenant, vous allez aligner ces reads nettoyés sur le génome de référence à l’aide de bwa [4] et samtools [5].

Combien de reads ne sont pas mappés ?

reads ne sont pas mappés.

Croisement de données

Calculez le nombre de reads qui chevauchent avec au moins 50% de leur longueur le gène trmNF grâce à l’outil bedtools [6]:

reads chevauchent le gène d’intérêt.

Visualisation

Utilisez IGV [7] sous sa version en ligne pour visualiser les alignements sur le gène. Faites une capture d’écran du gène entier.

References

1. toolkit NS. NCBI SRA toolkit. NCBI, GitHub repository. 2019.
2. Andrews S. FastQC a quality control tool for high throughput sequence data. http://www.bioinformatics.babraham.ac.uk/projects/fastqc/. http://www.bioinformatics.babraham.ac.uk/projects/fastqc/.
3. Zhou Y, Chen Y, Chen S, Gu J. Fastp: An ultra-fast all-in-one FASTQ preprocessor. Bioinformatics. 2018;34:i884–90. doi:10.1093/bioinformatics/bty560.
4. Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. arXiv preprint arXiv:13033997. 2013.
5. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9.
6. Quinlan AR, Hall IM. BEDTools: A flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26:841–2.
7. Thorvaldsdóttir H, Robinson JT, Mesirov JP. Integrative genomics viewer (IGV): High-performance genomics data visualization and exploration. Briefings in bioinformatics. 2013;14:178–92.
 

A work by Migale Bioinformatics Facility

https://migale.inrae.fr

Our two affiliations to cite us:

Université Paris-Saclay, INRAE, MaIAGE, 78350, Jouy-en-Josas, France

Université Paris-Saclay, INRAE, BioinfOmics, MIGALE bioinformatics facility, 78350, Jouy-en-Josas, France